On customer flows in Jackson queuing networks 



Sen Tan and Aihua Xia* 
Department of Mathematics and Statistics 
The University of Melbourne 
Parkville, VIC 3052 
Austraha 

19 January, 2010 



Abstract 

Melamed's theorem states that for a Jackson queuing network, the equihbrium flow 
along a link follows Poisson distribution if and only if no customers can travel along 
the link more than once. Barbour & Brown (1996) considered the Poisson approximate 
version of Melamed's theorem by allowing the customers a small probability p of trav- 
elling along the link more than once. In this paper, we prove that the customer flow 
process is a Poisson cluster process and then establish a general approximate version 
of Melamed's theorem accommodating all possible cases of < p < 1. 

Key words and phrases. Jackson queuing network. Palm distribution, Poisson cluster 
process, over-disperson. Stein's method, negative binomial. 
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1 Introduction 



We consider a Jackson queuing network with J queues and the following specifications [see 
Barbour & Brown (1996) for more details]. First, we assume that customers can move from 
one queue to another as well as can enter and leave from any queue. We assume that the 
exogenous arrival processes are independent Poisson processes with rates z^j, 1 < j < J. 
Service requirements are assumed to be exponential random variables with parameter 1 
and when there are m customers in queue j, the service effort for queue j is (f)j{m), where 
(f)j{0) = 0, 0j(l) > and (pjim) is a non- decreasing function of m. Second, we define the 
switching process as follows. Let A^^ be the probability that an individual moves from queue 
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i to queue j, /Xj be the exit probability from queue i and it is natural to assume 

J 

^ + /ii = 1, 1 < i < J. 

3=1 

Without loss of generality, we may assume that the network is irreducible in the sense that 
all customers can access any queue with a positive probability. Set aj as the total rate of 
arriving customers (including both exogenous and endogenous arrivals) to queue j, then the 
rates {c(j} satisfy the equations 

J 

i=l 

and they are the unique solution of the equations with aj > for all j. 

For convenience, we define state as the outside of the network, that is, the point of 
arrival and departure of an individual into and from the system. We write S := {{j,k) : 
< i, < </} as the set of all possible direct links and use H-'^ to record the transitions 
of individuals moving from queue j to queue k, then S = {S-''^, < j,k < J} gives a 
full account of customer flows in the network, where departures are transitions to and 
arrivals are transitions from 0. If pjk is the rate of equilibrium flow along the link (j, k), then 
Pjk = ctj^jk and the mean measure of S is 

\{ds, (j, k)) = pjkds, sen, (j, k) e S. 

Our interest is on the customer flows along the links in C G S for the time interval [0, t], 
so we set the carrier space as Tct = [0, t]x C and use Hc,t to stand for the transitions along 
the links in C for the period [0,t]. Then the mean measure of Se t is 

><c,t{(^^^ 0'' ^)) = Pjkds, 0<s<t, {j, k) e C. 

Melamed's theorem states that H,^ ^ is a Poisson process if and only if no customers travel 
along the links in C more than once [Melamed (1979) proved the theorem when (f)j{m) = Cj 
for m > 1 and 1 < j < J, and the general case was completed by Walrand & Varaiya (1981)]. 
Barbour & Brown (1996) considered the Poisson approximate version of Melamed's theorem 
by allowing the customers a small probability of traveling along the links more than once. 
For convenience, we call the probability of customers traveling along the links in C more 
than once as the loop probability. The bounds for the errors of Poisson process approximation 
are sharpened by Brown, Weinberg & Xia (2000) and Brown, Fackrell & Xia (2005) and it 
is concluded in these studies that the accuracy of Poisson approximation depends on how 
small the loop probability is. 

In section 2 of the paper, we use the Palm theory, Barbour-Brown Lemma [Barbour 
& Brown (1996)] and infinite divisibility of point processes to prove that is a Poisson 
cluster process [Daley & Vere- Jones (1988), p. 243]. The characterization involves a few 
quantities which are generally intractable, so in section 3, we prove that H is over-dispersed 
[see Brown, Hamza & Xia (1998)], i.e., its variance is greater than its mean, and conclude 
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that suitable approximations should be those with the same property, such as the compound 
Poisson or negative binomial distributions. We then establish a general approximate ver- 
sion of Melamed's theorem for the total number of customers traveling along the links in 
C for the period [0,t] based on a suitably chosen negative binomial distribution. The ap- 
proximation error is measured in terms of the total variation distance and the error bound 
is small when the loop probability is small [cf. the Poisson approximation error bound in 
Barbour & Brown (1996)] and/or t is large with order [cf. Berry-Esseen bound for normal 
approximation, Petrov (1995)]. 



2 A characterization of the customer flow process 

The customer flows are directly linked to the changes of the states of the queue lengths, so 
we define Ni{t) as the number of customers, including those in service, at queue i and time 
t, 1 < i < J . Then {N(t) := {Ni{t), . . . , Nj{t)) : t G IR} is a pure Markov jump process on 
state space {0, 1, 2, ... j"' with the following transition rates n = (rii, . . . , rtj): 

qn,n+ej = Vj, qn^n-ej+e^ = (j)j{nj)Xjk, qn,n-e, = fij(f)j{nj), 1 <j,k < J, 

where ej is the jth coordinate vector in {0, 1,2,... . The stationary Markov queue length 
process has a unique stationary distribution: for each t > 0, Nj(t), j = 1, . . . , J are indepen- 
dent with 



F{Nj{t) = k) 



Let Xi be the ith queue visited by a given customer, then {Xi : z = 1, 2, . . . } is called 
the forward customer chain and it is a homogeneous finite Markov chain with transition 
probabilities 



Poo = 0, Pok = J 
Pjo = /Wj, Pjk = \k, j,k = l,...,J. 



The backward customer chain X* is the forward customer chain for the time-reversed process 
of {N(t) := (iVi(t), . . . ,iVj(t)) : t e IR} [Barbour & Brown (1996), p. 475] and it can be 
viewed as the time- reversal of the forward customer chain {Xi : i = 1,2, . . .} with transition 
probabilities 



pSo = 0' Poi 



* _ * _ %Ajfc ■ L _ 1 T 

PkO ~ ' Pkj ~ ) ]) — ^) ■ ■ ■ ) 'J- 

ak ak 

We will use the Palm distributions to characterize the distribution of Sc^^, prove its prop- 
erties and establish a general approximate version of Melamed's theorem for 'B^ctiJ'c.t)- For 
the point process H with locally finite mean measure \{ds, (j, k)) = pjkds, s G IR, (j, k) G S, 
we may consider it as a random measure on the metric space IR x 5 equipped with the metric 

diiui, (ji, j2)), (m2, iki,k2))) = \ui - M2|loij2)^(fci,fc2)> Ui,U2 G R and (ji, j2), iki,k2) G S, 
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so that we can define the Palm distribution at a ElRxS as the distribution of S conditional 
on the presence of a point at a, that is, 

pa . _ I [^e] V — a^Yct>^- almost surely, 
A(cta) 

see Kallenberg (1983), p. 83 for more details. A process 3" is called the Palm process of S 
at a if its distribution is P". In applications, it is often more convenient to work with the 
reduced Palm process S° — 6a [Kallenberg (1983), p. 84], where Sa is the Dirac measure at a. 

The Palm distributions are closely related to the size-biasing in sampling contexts 
[Cochran (1977)]. More precisely, if X is a non-negative integer- valued random variable, 
one may consider it as a point process with the carrier space having only one point so its 
Palm distribution becomes 

^ ^ EX 

However, this is exactly the definition of the size biased distribution of X [see Goldstein & 
Xia (2006)]. 



Lemma 2.1 [Barbour & Brown (1996)] For the open queuing network, the reduced Palm 
distribution for the network given a transition at link (j, k) at time is the same as that 
for the original network, save that the network on (0, oo) behaves as if there were an extra 
individual at queue k at time and the network on (— oo, 0) behaves as if there were an extra 
individual in queue j at time 0. 



For two random elements rji and 772 having the same distribution, we write for brevity 

d 



Lemma 2.2 For each [j, k) G S , there is a point process on]R x S independent of 

H such that 

^(o,(j,fc)) + 3 A 3(o,(i,fc))_ 



Proof. The proof is adapted from Barbour & Brown (1996), p. 480. By Lemma l2.ll the 
reduced Palm process H'-'^'^-'''^-'^ — 5(o,(j,fc)) has the same distribution as that of H except that 
the network on (0, 00) behaves as if there were an extra individual at queue k at time 
and the network on (—00, 0) behaves as if there were an extra individual in queue j at time 
0. Let and X^'^^ be the routes taken by the extra individual on (— oo,0) and (0, 00) 
respectively. Whenever the extra customer is at queue i together with other m customers, 
we use independently sampled exponential service requirements with instantaneous service 
rate 0i(m + 1) — (f)i{m). Noting that this construction ensures that the extra customer uses 
the "spare" service effort and never "interferes" with the flow of the main traffic, one can see 
that its transitions are independent of 3. The same procedure applies to the construction of 
the backward route. Let be the transitions taken by the extra customer on (—00, 0)U 

(0, 00) plus the Dirac measure 5(o,(j,fc))5 then is independent of S and the conclusion 

of the lemma follows from the construction. ■ 
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Let 6s, s G R, denote the shift operator on R x 5 which translates each point in R x 5 
by s to the left, i.e. Os{{u, (j, k))) = [u — s, (j, k)) and we use to stand for a copy of 

^(o,0»)o^„sGR. 

From now on, we focus on the point process Sc,*. With metric d, Tct is a Polish space 
and we use B (Tct) to stand for the Borel a-algebra in Tc,t- Let Hc^t denote the class of all 
configurations (finite nonnegative integer-valued measures) on Tc^t with Hct the a-algebra 
in Hc^t generated by the sets 

{e e Hc,t : as) =^}, Z+ := {0, 1, 2, . . . }, 5 G /5 {Tc,t) , 
see Kallenberg (1983), p. 12. 

Theorem 2.3 Let {r/^, i > 0} be independent and identically distributed random measures 
on Tc.t having the distribution 



P [rio i^c,t) > 1] = 
where 




Let M be a Poisson random variable with mean 9c,t o,nd independent of {rj^, i > 0} , then 

M 
i=l 

Proof. By Lemma [212] and Theorem 11.2 of [Kallenberg (1983)], we can conclude that H^t is 
infinitely divisible, hence we obtain from Lemma 6.6 and Theorem 6.1 of [Kallenberg (1983)] 
that Se t is a Poisson cluster process, that is, 

M 
i=l 

where r/j, i > are independent and identically distributed random measures on Tc,t such 
that P (*7o(rc,t) > 1) = 1, M is a Poisson random variable with mean 9c,t and independent 
of {rfi, i > 1}. The direct verification ensures that the Palm process of ^^^i ry^ at a e Tct is 
Z^i=i Vi + ^0 ) where tJq is the Palm process of tiq at a and is independent of {M, r]^, i > 1}. 
This in turn implies that ^^'^^^^''^^ = rjl^''^^'''^\ 

Let fi{ds, (j, k)) denote the mean measure of the point process rj^, then some elementary 
computation ensures that the mean measure of "^fli »7j is 9c,tfJ'{ds, (j, k)) for (j, k) E C and 
< s < t. On the other hand, the mean measure of Sc,t is \c,t{ds, (j, k)) = pjkds, (j, k) E C 
and s G [0,t], so we obtain 

^i{ds, (j, k)) = §^ds, (j. A;) G C, s G [0, t]. (2.3) 
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The representation (12.11) is because of the fact that P {r)Q{rc,t) > 1) = 1 and 

Jtc, rio (Tct) (^.^^ Jo ^^''^^''^^ (Tct) ^c,t 

In particular, if we take A = Hct, then the left hand side becomes 1, so ( 12 ■2p follows. ■ 

Despite the fact that 6c^t is specified by (12.21) . since the Palm process ^^^'^^'''^^ is generally 
intractable, it is virtually impossible to express 6c,t explicitly in terms of the specifications 
of the Jackson queuing network. On the other hand, the relationship (12. 3p yields 

Er7o(rc,) = ^t. 

fc,t 

The following proposition tells us the range of values that Erj^^Tct) and 6c\t niay take. To 
this end, we define 

6c;(j, A;) = E^(°'(^'^»(PxC)-l and ec= Yl —^c{j,k). (2.4) 

(i,fc)ec 

In other words, ec{j, k) is the average number of visits in C by the extra customer crossing 
the link (j, k) and ec is the weighted average number of visits by an extra customer crossing 
links in C. 

Proposition 2.4 We have 

l<Er7o(rc,t)<l + ec (2.5) 

and 

t<ec,t<pct. (2.6) 



Proof. The first inequality of (12. 5p follows immediately from the fact that 
'P{r}Q{Vc,t) > 1) = 1- For the second inequality of (12. 5p . noting that the mean measure of ryg 
is 



we have 



2 

PC ■ 



c,t 



= [E77o(rc,)]' < E [r7o(rc,)^] 



U,k)&C • 

< Y.(^+'cij,k))l^t 

ij,k)ec 

PC ^ ^ ^{j,k)ec^c{j,k)pjk ^ 



We divide both sides by -^t to get 

Er7o(rc,t) = ^t<l + ec. 
^c,t 

Finally, fl2.6p is an immediate consequence of (12 .Sp and the equation Oc^t = Brf%ct) ^' " 

Remark 2.5 If the loop probability in C is 0, then ec = and lE77o(rc t) = 1, so there is 
only one customer on Tc,t- This customer is crossing the link (j, k) with probability ^ at a 
time uniformly distributed on [0,t]. 

3 A discrete central limit theorem for the customer 
flow process 

A random variable is said to be over- dispersed (resp. under-dispersed) if its variance to mean 
ratio is greater (resp. less) than one. A random measure x on a Polish space is said to be 
over- dispersed (resp. under-dispersed) if x{B) is over-dispersed (resp. under-dispersed) for all 
bounded Borel subset B of the Polish space. It is concluded in Brown, Hamza & Xia (1998) 
that point processes which arise from Markov chains which are time-reversible, have finitely 
many states and are irreducible are always over-dispersed. As our process H is virtually a 
multivariate version of point processes studied in Brown, Hamza & Xia (1998), the following 
property can be viewed as a natural extension of the study in Brown, Hamza & Xia (1998). 

Proposition 3.1 The point process 3 is over- dispersed. 

Proof. The space (IR x S, d) is a Polish space and for each bounded Borel subset B of 
R X iS, it follows from the definition of the Palm processes [see Kallenberg (1983), p. 84, 
equation (10.4)] that 

1E[H(5)]^ = E I H(5)H(cia)=E I H°(5)A(da) > E I (H(5) + 1) A(cia), 
JB Jb Jb 

that is, 

Var [H(5)] > EH(5), (3.1) 

completing the proof. ■ 

The inequality in (13. ip is generally strict except that the loop probability is 0, i.e. H is a 
Poisson process. Hence, suitable approximate models for the distribution of Sc,i := Sc,t(Xc,t) 
are necessarily over-dispersed. One potential candidate for approximating the distribution of 
Set is the compound Poisson distribution. However, as it is virtually impossible to extract 
the distribution of ^^'^'^•^'''^\Tc\t) for < s < t, we face the same difficulty to specify and 
estimate the approximate distribution if a general compound Poisson is used. On the other 
hand, as a special family of the compound Poisson distributions [Johnson, Kemp & Kotz 
(2005), pp. 212-213 and p. 346], the negative binomial distribution has been well documented 
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as a natural model for many over-dispersed random phenomena [see Bliss & Fisher (1953), 
Wang & Xia (2008) and Xia & Zhang (2009)]. The negative binomial distribution NB(r, g), 
r > 0, < g < 1, is defined as 

T(r + i) ^.^ .j . ^ 
i [rjtl 

The advantage of using negative binomial approximation is that it suffices to estimate the 
mean and variance of the approximating distribution, like what we often do in applying the 
central limit theorem based on the normal approximation. 

We will use the total variation distance between the distributions of nonnegative integer- 
valued random variables Yi and I2 

drviYuY^) := sup |P(ri e A) - P(F2 G A)\ 

to measure the approximation errors in negative binomial approximation. 

The discrete central limit theorem is valid under the assumption that the loop probability 
is less than 1. More precisely, let wc{jk) be the probability that a link (j, k) crossing customer 
crosses the links in C only once, i.e., the only time that the customer crosses the links in C 
is the one the customer is crossing. Define 

wc= wc{3k)pjk/pc, 

the weighted probability of customers crossing links in C only once. Clearly, we have wc{jk) > 
fXk, so 

{j,k)£C 

The following lemma plays a crucial rule for the estimation of the negative binomial approx- 
imation error. 

Lemma 3.2 drv {^c,t, ^c,t + 1) < ^i=^- 

Proof. We prove the claim by a coupling based on the "priority principle" [cf. the proof 
of Lemma 12.2] . We define a customer as a single crossing (sc for brevity) customer if the 
customer crosses links in C only once, otherwise, the customer is labeled as multiple crossing, 
or mc for short. We "manage" the network by regrouping the customers at each queue into 
sc customers and mc customers. Whenever there are m2 mc customers together with mi sc 
customers at queue j, we use independently sampled exponential service requirements with 
instantaneous service rate 0j (mi -1-777.2) ~0j('^i) for all of the mc customers while the service 
for the sc customers is carried out with instantaneous service rate (j)j{mi), that is, as if no mc 
customers present at the queue. Since the sc customers take priority over the mc customers 
and the 777c customers use the "spare" service effort and never interrupt the traffic fiow of the 
sc ones, we can see that its transitions are independent of the transitions of the sc customers. 
Let (resp. Z^^) denote the transitions of sc (resp. mc) customers moving from queue j 
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to queue k in the period [0,t], then Zi := {ZC^, [j, k) G C} and Z2 := {^2^, (j, k) G C} are 
independent and 

3c,t = Zi + Z2. 

By Melamed's theorem, the point process Zi is a Poisson process with mean measure 

Azi(c^s, (j, k)) = wcijk)pjkds, (j, /c) G C, 0< s <t, 
so Zi(rc7() follows Poisson distribution with mean wcPct and 

1 



drv (Sc,t, Sc,t + 1) < drv C^iiT c,t), M^a) + 1) < 



^2ewcPct'' 



where the last inequality is because of the fact that the distribution of Poisson is unimodal 
and Proposition A. 2. 7 of [Barbour, Hoist & Janson (1992), p. 262]. ■ 



To state the discrete central limit theorem, we set 



k) = E[^(°'(^''=))(1R X C)(|(°'(^''=»(R x C) - 1)] and ac = —""cU, k). 

{j,k)&C 

That is, (Jciii k) is the second factorial moment of the number of visits in C by the extra 
customer crossing the link (j, k) and ac is the weighted average of the second factorial 
moments of the number of visits by an extra customer crossing links in C [cf. f l2.4p ]. 

Theorem 3.3 Let 

{pctf pet 



Var(Sc,t) - Pet ' Var(Sc,t) ' 

then 

dTy(Sc,t,NB(r,g)) < ] {2(Var(Sc,t) - Pct)^ 

{pctY\/2ewcpct 

+pct(Sc,43] - Pcf^cA^] - 2pctiVaTiEc,t) - Pet))} (3.2) 

< ^=i==(26^ + ac), (3.3) 
^/2ewcPct 

where Sc.tln] stands for the nth factorial moment of^ct defined as 

Sc^N = ]E[Sc,i(Hc,t - 1) . . . {^c,t - n + 1)]. 

Remark 3.4 The parameters of the approximating negative binomial distribution are cho- 
sen so that it matches the mean and variance of '^c,t- 

Remark 3.5 If the loop probability in C is 0, then the negative binomial is reduced to 
Poisson distribution and the upper bound in Theorem 13.31 becomes 0. This implies half of 
Melamed's theorem (1979). 
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Remark 3.6 If the loop probability is between and 1, then both ec and ac are finite, so the 
negative binomial approximation error bound is of order 0(1/ \/t). Furthermore, if the loop 
probability is small, then both ec and ac are small, so the negative binomial approximation 
to the distribution of is even faster. 



Proof of Theorem 13.31 The essence of Stein's method is to find a generator which charac- 
terizes the approximating distribution, establish a Stein identity to transform the problem 
of estimating the approximation errors into the study of the structure of the object under 
investigation. In the context of negative binomial approximation, let a = r(l — g), b = I —q, 
then a generator which characterizes NB(r, g) is defined as 

Bg{i) = (a + bi)g{i + 1) — ig{i), i G Z+, 

for all bounded functions g on Z+ [see Brown & Phillips (1999) and Brown & Xia (2001)]. 
The Stein identity is naturally established as 

Bg{^ = fii) - TiU) (3.4) 

for / G J-" := {/ : Z+ — )• [0,1]}, where rc{f) = J2°^q f{i)'n'i- It was shown in Brown & 
Xia (2001) that, for each f & J^, the solution gf to the Stein equation (13. 4p satisfies 

l|A^?/|| < -, (3.5) 
where Agf{-) = gf{- + 1) — gf{-). The Stein identity (13. 4p ensures that 

sup |E/(Sc,t) - vr(/)| = sup \EBgf{Ec,t)\ , 
feT far 

hence, it suffices to estimate WjEgf^Ect) for all f E J^. For convenience, we drop / from the 
subindex of gj. By Lemma 12.21 we can take a point process $}ct^'^^^ on independent of 
Sc^f such that 

Therefore, if we write ^ci^'''^\Tc,t) = 1 + then 

E%(Hc,t) = nia + bEc,t)9{^c,t + l)-^cM^c,t)] 

= alEgiEc,t + l) + b f 9{^c,t + 2 + i^''^^'^^^)p,uds 

J2 f g{Ec,t + l + &'^^''^^)Pjkds. (3.6) 



(i,fc)GC • 



Let 



a + (6-1) J2 Pjkt = 0, (3.7) 
(i,fc)ec 
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and Hc,t = Hc,t + 1, then it follows from fl3.6p that 



{j,k)ec 
(i,fc)ec 



Pifc 



r=0 



bAg Ec,t > PjfcC^s 



(3. 



Now, set 



^imec lo n^'^^'^'^^P^kds _ Var (Ec,) - Pet 



E(, .)ec lo Ee(^'(^>))p,.rfs + pet Var {Ec,t) 



j{j,k)£C Jo 

where the last equality is due to the following observation: 



(3.9) 



c,t 



V E / {Ec,t + l + &^'''^^)pjkds 



■''0 



(i,fc)ec • 



and so 



{j,k)&C 

We then obtain from fl3.8D that 



E T [\^''^''''^^Pjkds = V^T{Ec,t)-pct. 



(3.10) 



MBg{E 



c,t 



(i,fe)ec 

E E 

(i,fc)6C 





'g(s,(3,fc))_l - 




E 


r=0 








E 


r=0 



1=0 
r-l 



Pjkds 



b]EA'g {Ec,t + r) - (1 - 6) 5^ EA^^^ {Ec,t + I 



Pjkds, 
(3.11) 



where the last equation is due to the fact that ^(^'^-J'*^)) is independent of E^t- On the other 
hand, using (13. 5p . we have 



lEA^giEct + 1) < 2\\Ag\\dTvi^c,u ^c,t + 1) < 



2dTv{'^c,t, '^c,t + 1) 
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so it follows from (13. lip that 



|E%i 



-c,t) 



< ^^""^'f^'' ^ E / [26Ee(^'(^>)) + (1 - 6)Ee(^'(^''=)) - 1)] p,^^^^^ 

(3.12) 

Using the Palm distributions of Set together with (I3.10p . we get 

Sc,i[3] = e/ {Ec,t-l){^c,t-2)Ec,t{da) 

= E r E [i^c,t + &^^^^'^^) (He, + - 1)] 

(i,fc)ec 

= pctHc,42] + 2pct(Var(Hc,0 - Pet) + V /" E^^^'^^'^^" (C^^'^^''^" - l) p,,ds. 



(i,fc)6C 

This in turn ensures 



f Ee^^'^^''^" (e^^'^^''^^ - 1) P.kds = EcA3]-pctEcA2] -2pct{VaT{Ec,t)-pct). (3.13) 



Consequently, combining (1330D . (137[3D with (1332|) gives ([32]). 
Finally, by the definitions of ec and (Xc, we have 



and 



J2 I 'E^^'''^^''''>\i'^''^^'''^'> -l)p,uds<acpct. 
{j,k)ec 

Therefore, ([SJ]) follows from ([321), (EHOD and (Km . m 
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